Understanding exoplanet atmospheres

# Import required libraries
import numpy as np
import pandas as pd

When a planet passes in front of its star, the tiny dip in brightness leaves a fingerprint. That transit depth, measured across wavelengths, is what this project digs into. Each observation carries noise, so the real question is bigger than what does the spectrum look like? — it’s what physical properties of the planet actually drive that signal, and how much can we trust our reading of them? That’s the lens I use throughout: interpretability and sensitivity.

# Load the CSV file into a Pandas DataFrame
df = pd.read_csv("small_astro.csv")

df_backup = df.copy()
df_backup.head(10)

# Columns of interest for planetary and chemical properties and transit depth
feature_columns = df.columns[1:8]
transit_depth_column = 'x1'  # pick the first wavelength

# Calculate mean and variance for features and transit depth
feature_mean = df[feature_columns].mean()
feature_variance = df[feature_columns].var()
transit_depth_mean = df[transit_depth_column].mean()
transit_depth_variance = df[transit_depth_column].var()

feature_mean, feature_variance, transit_depth_mean, transit_depth_variance

What the data looks like

The dataset carries a handful of planetary and chemical features — radius, temperature, and the log concentrations of H₂O, CO₂, CO, CH₄, and NH₃ — alongside the transit depth at a chosen wavelength. The spreads are informative. Planet radius hovers around 0.70 with modest variance, while temperature swings over a much wider range. The gas concentrations are all low, as you’d expect in log space, and the transit depths cluster near 0.01. Whatever drives that signal, it isn’t a single obvious spread — which is exactly why sensitivity analysis matters.

# Import required libraries for modeling and SHAP values
import xgboost as xgb
import shap

# Prepare the feature matrix (X) and the target vector (y)
X = df.iloc[:, 1:8]
y = df['x1']

# Train an XGBoost model
model = xgb.XGBRegressor(objective='reg:squarederror')
model.fit(X, y)

# Initialize the SHAP explainer and compute contributions
explainer = shap.Explainer(model)
shap_values = explainer(X)

Attributing predictions with SHAP

SHAP comes from cooperative game theory. Think of each feature as a player in a game and the prediction as the payoff: a feature’s Shapley value is the average change in the payoff when that player joins the coalition. I train an XGBoost regressor on the features and compute SHAP values for every observation, which tells me not just what the model thinks is important but how each feature pushes a given prediction up or down.

# Import the Random Forest Regressor
from sklearn.ensemble import RandomForestRegressor

# Train a Random Forest model
rf_model = RandomForestRegressor(n_estimators=100, random_state=42)
rf_model.fit(X, y)

# Extract feature importances
feature_importances = rf_model.feature_importances_

# Look at the importance ranking
importance_df = pd.DataFrame({
    'Feature': feature_columns,
    'Importance': feature_importances
}).sort_values(by='Importance', ascending=False)
importance_df

Random forest importance

As a second, more direct read on importance, I fit a random forest regressor and inspect its built-in feature importances. The picture is consistent: planet radius is the dominant driver of transit depth, with the CO₂ and CO concentrations carrying most of the remaining signal. Temperature, by comparison, barely moves the needle.

Cross-checking with a sensitivity metric

I also bring in a model-agnostic sensitivity measure called proportional marginal effects (PME). It approaches the question from a different direction than SHAP or permutation importance but lands on a similar conclusion: planet radius stands out as the pivotal feature, with the gas concentrations trailing and temperature looking marginal. When several tools tell you the same thing, that’s a good sign it’s real rather than an artifact of one method.

What this means for uncertainty and parsimony

The practical payoff is twofold. First, uncertainty: a feature with high sensitivity means measurement error in it propagates directly into the prediction. If planet radius is the dominant driver and it’s hard to measure precisely, that uncertainty belongs in the final forecast. Second, parsimony: when a model has a clear hierarchy of importance, you can drop the marginal features and simplify without losing much. That reduces overfitting risk in small data and tells observers where to spend their measurement effort. The trade-off, of course, is that simplifying always costs some predictive accuracy — the goal is to do it knowingly.

A weighted view with PCA

I wanted a low-dimensional view of the dataset that still respects what matters, so I run a weighted PCA: standardize the features, scale each one by the square root of its importance, then project onto the first two principal components. The weighting nudges the projection toward the structure driven by the important features, and the explained-variance share tells me how much of the original structure survives in two dimensions.

from sklearn.decomposition import PCA
from sklearn.preprocessing import StandardScaler
import numpy as np

# Standardize the feature matrix
X_scaled = StandardScaler().fit_transform(X)

# Weight each feature by its importance (sqrt) for weighted PCA
X_weighted = X_scaled * np.sqrt(feature_importances)

# Perform PCA on the weighted data
pca = PCA(n_components=2)
X_pca = pca.fit_transform(X_weighted)

# Variance explained by the first two principal components
pca.explained_variance_ratio_

Reference

Changeat, Q., & Yip, K. H. (2023). ESA-Ariel Data Challenge NeurIPS 2022: Introduction to exo-atmospheric studies and presentation of the Atmospheric Big Challenge (ABC) Database. arXiv preprint arXiv:2206.14633.

Herin, M., Il Idrissi, M., Chabridon, V., & Iooss, B. (2022). Proportional marginal effects for global sensitivity analysis. arXiv preprint arXiv:2210.13065.